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ABSTRACT 


An experimental and analytical investigation was conducted to develop a calibrated 
torsional model of a three-cylinder, two-stroke diesel engine. A Detroit Diesel 3-53 engine 
was instrumented for time resolved measurement of cylinder firing pressures and high 
resolution near instantaneous shaft speed using a 720 and a 3,600 count per revolution optical 
encoder. Data were taken for three speeds and three torques for a total of nine conditions. A 
SIX Foci feedoti torsional vibration model of the crankshaft, connecting rods, and 
pistons was developed. The nonlinear inertias, due to the reciprocating pistons, were included 
along with linear stiffness and damping. The equations of motion were numerically integrated 
over a cycle to obtain predicted response. The predicted response was compared to the 


measured response at the free end of the crankshaft. 
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I. INTRODUCTION 


A. MOTIVATION 

Diesel engines are widely used in both military and civilian applications. In the 
Navy, they are utilized for main propulsion and electric power generation. The reliable 
detection of faults both non-obtrusively and inexpensively could shift the focus of 
maintenance from regularly scheduled to condition based. This maintenance strategy has 
the potential to save costly, unneeded maintenance and can increase availability of diesel 
engines. 

Individual cylinder pressures of internal combustion engines are excellent 
predictors of many types of engine faults. Cylinder pressures are directly related to the 
torque output capability of an engine. Therefore, the contribution of torque from an 
individual cylinder is directly related to the maximum pressure that cylinder develops. 
Cylinders experiencing loss of pressure can be caused by faulty fuel injectors, worn 
cylinder lining, faulty piston rings, leaking valves, and blown or leaking cylinder head 
gaskets. 

There is a problem with measuring cylinder firing pressures in the field. While high 
quality pressure transducers can accurately measure cylinder pressures, they have a limited 
lifetime in the combustion chamber environment and are expensive. Many alternative 
indirect methods of cylinder pressure monitoring such as the use of strain bolts are being 


explored. 


Another indirect method is to measure the variations in shaft speed. With this 
information and equations of motion describing the dynamics of the engine, the firing 
pressures are solved for numerically. In order to accurately predict cylinder firing 
pressures with shaft speed variation, an accurate torsional engine model representative of 


the physical system is required. 


B. STATE OF THE ART 

There have been numerous research efforts focused on extracting cylinder pressure 
information from crankshaft speed fluctuations. It 1s well known that these speed 
fluctuations are a direct result of the cyclic nature of the crankshaft torque waveform. 
This torque waveform is a direct result from gas pressure forces, due to combustion, and 
inertia forces, due to internal reciprocating masses. The challenge has been taking these 
measured speed fluctuations and passing them through a robust engine model to predict 
cylinder pressures. This literature review focused in areas of non-obtrusive cylinder 
pressure and torque prediction with emphasis on torsional vibration modeling. 

Wilson [Ref. 1] developed a torsional vibration model of a crankshaft with 
reciprocating components. With his methods, the individual inertia, stiffness, and damping 
between cylinder crankthrows and the rest of the shaft are readily determined. His 
methods are currently used by industry to develop torsional vibration models of 
crankshafts. Kabele [Ref. 2] developed a torsional vibration model of a diesel crankshaft 
which included piston and ring friction, connecting rod plane motion, hysteretic damping 


in the crankshaft, and cylinder pressure as a forcing function. The model was derived 


from Newtonian Mechanics and the governing differential equations were non-linear 
making a numerical solution necessary. Once his model was developed, he compared 
measured and model predicted vibration amplitudes for various harmonic orders and 
concluded that correlation between the two was excellent. The measured vibration 
amplitudes were from a V-8 diesel engine. 

Citron et al, [Ref 3] did a computer simulation and produced cylinder pressure 
torque waveforms by using an elastic model of an engine, drivetrain system, and measured 
speed fluctuations. These speed fluctuation data were passed back through the model to 
determine the fluctuation waveform of both the total engine torque being developed and 
the cylinder pressure waveform which gave rise to it. 

Rizzoni [Ref. 4] proposed passing crankshaft speed fluctuations through an 
equivalent electronic circuit representing engine dynamics to predict cylinder misfires. He 
conducted experimental work on production spark ignited vehicles and confirmed that it is 
possible to apply this method with minimal hardware and computational overhead. 

Sobel et al, [Ref. 5] described a technique for measuring instantaneous crankshaft 
torque utilizing a non-contact ferromagnetic material to measure torsional stress in the 
crankshaft. The instantaneous torque of an internal combustion, four-cylinder engine was 
measured at the flywheel and compared with pressure signals during bench tests. 
Excellent correlation between this non-contact measured torque and values of the mean 
effective pressure for each cylinder were obtained. 

Brown and Neill [Ref 6] presented a non-contact method for determining 


simultaneously the pressure in each cylinder by employing a pattern recognition technique 


to compare crankshaft speed fluctuations to reference patterns in a knowledge base. The 
non-contact method utilized was an interval timer and a magnetic sensor which timed the 
flywheel gear teeth as they passed. The experiment was carried out on a Detroit Diesel 
6V-92TA engine. The experimental results show that the method evaluates the cylinder 
pressures with an RMS error of less than six percent. But this requires extensive data 
bases for each condition. 

Bell [Ref. 7] concluded, utilizing the same Detroit Diesel 3-53 engine used in this 
research, that the variation of shaft speed held information which could be used to predict 
cylinder firing pressures. He also recommended the development of a torsional vibration 
model of the engine and flywheel system along with a more rigid optical encoder mounting 


for more accurate time resolved measurement. 


e METHODOLOGY AND OBJECTIVES 

The engine used for this research is a Detroit Diesel Series 53 engine, model 5033- 
SOOIN. it is a three-cylinder, two-stroke engine. During one revolution, the crankshaft 
speed fluctuates as a result of cylinder pressure variations. These crankshaft speed 
fluctuations can be measured using an optical encoder. A torsional vibration model with 
nonlinear inertias, due to the reciprocating pistons, and linear stiffness and damping will be 
developed to represent the dynamics of the engine. The torsional vibration model will be 
calibrated by ensuring, as closely as possible, predicted response matches the measured 
response. Data from these measured speed fluctuations can then later be passed through 


the calibrated torsional model to determine cylinder firing pressures 


The objectives of this thesis are to: 
1. Instrument a diesel engine with a rigid optical encoder mount such that reliable time 
resolved measurements of crankshaft speed relative to the engine block can be obtained. 
2. Develop a torsional vibration model of the crankshaft including components connected 
to it (i.e. connecting rods, pistons, flywheels, and external loads). 
3. Operate the diesel engine at several speeds and applied torques to establish baseline 
data. 


4. Compare the model predicted response to the experimental measured response. 


D. ORGANIZATION 

Chapter II describes the development of the crankshaft torsional vibration model. 
Equations of motion are derived for a six degree-of-freedom torsional model and 
numerical solution procedure is presented. 

Chapter II] describes the engine, dynomometer, instrumentation, and data 
acquisition procedures. 

Chapter IV presents results of data acquisition, the measured response, and the 
model predicted response. A comparison of measured to predicted response is made 
based on the results. 


Chapter V contains summary, conclusions, and recommendations for future work. 


se 





Il. TORSIONAL VIBRATION MODEL 


A. PHYSICAL SYSTEM 

The engine used for this research is representative of diesel engines. It is a three- 
cylinder, two-stroke engine. Engine operating conditions and design specifications are 
discussed in more detail in Section III. A. 

Figure 2.1 illustrates the engine layout with the crankshaft, flywheel, pistons and 
connecting rods. The front of the engine is to the left. The cylinders are numbered from 
front (1) to back (3). Connected to the power take-off shaft on the mght is the load, 
which in this case is supplied by a water brake dynamometer. The firing order 1s 1-3-2 and 


the crankshaft rotates clockwise as seen from the front. 


INJECTORS 





P ie. 
VALVE COVER 


: | “HEAD ASSEMBLY 


4a 
ENGINE BLOCK 


IES wa Take OX¢ 


® | : 
Front - i : } . e : 





DETRBIT DIESEL 53° SERIES 


Figure 2.1 Detroit Diesel 3-53 with crankshaft and major reciprocating 
components. The front of the engine is the left hand side. 
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The crankshaft crankthrows are the locations where the connecting rods are attached 
to the crankshaft. Crankthrows number one, two, and three are displaced 120 degrees. 
After the number one piston reaches top dead center (TDC), its highest point of travel, 
120 degrees of crankshaft rotation later the number three piston arrives at its TDC. 
Number two piston’s TDC is 240 degrees after number one. The power into the 
crankshaft is delivered at each crankthrow. 

The forces due to gas pressure in each cylinder can be taken as a force ‘-F’ on the 
piston and an equal and opposite force ‘F’ on the cylinder head. The force ‘-F’ is 
transmitted through the connecting rod to the crankthrow. The resultant force on the 
crankshaft multiplied by the eccentricity ‘R’ of the crankthrow creates the gas torque. If 
there were no friction, the gas torque would equal the shaft applied torque. Figure 2.2 


shows the relationship between gas pressure force and torque. 


Crankshaft 






—_——_— — — eee ee eee 


) Shaft Applied Torque 


Chamber 


where: F = cylinder pressure * area of the piston (Ibf) 
R = crankthrow eccentricity (in) 
8 = angular position of crankshaft relative to TDC (rad) 
L = length of connecting rod (in) 


Figure 2.2 Relationship between gas pressure forces and gas torque into the crankshaft. 
8 


Gas torque into the crankshaft as a result of cylinder gas pressures is given by the 


following relationship: 


T icy=PicyRsin( @ )A, where: Py-y; = Measured cylinder pressure (psia) 
Trey=Pr2eyiRsin( @ - 47/3) A, R = Crankshaft crankthrow eccentricity (in) 
T3cy=P3-yiRsin( @ -270/3)A, @ = Angular position of #1 crankthrow 
after TDC (rad) 
(2.1) A, = Cross-sectional area of the pistons (in’) 


The torque applied to the crankshaft is somewhat less than the gas torques due 
mainly to friction on the cylinder wall (rings). The major components considered in the 
development of a torsional vibration model are the crankshaft, pistons, connecting rods, 


and flywheel. This physical system to be modeled is shown in Figure 2.3. 


Number 2 Piston and 


| Connecting Rod. Crankthrow 
Crankthrow for Number 3 


for Number | 





Optical Encoder 
Attachment Location. 


Output End of 
Crankshaft where 
Flywheel and 
Load are Attached. 


RS 


The Four Joumal Bearing Surfaces 





Figure 2.3 Crankshaft with a single reciprocating component. 
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1. Model Development 

The reciprocating components of the engine are very complex. However, to 
facilitate calculation it is useful to represent the complex system by a reduced order 
system. For this research, six rigid masses connected by sections of shafting with torsional 
stiffness and damping shall be used. This reduced system should retain the dominant 
dynamic characteristics of the original physical system that are of interest. The six degree- 


of-freedom torsional vibration model that is used to represent the physical system is shown 


in Figure 2.4. Piston/Crank 
Number 3 Flywheel 


Piston/Crank 
Number 2 





co eo & 
rey Guy 
ee ee ee ee ee 
Ky ia, 
T, 
qa 
Cz 





Figure 2.4 Torsional vibration model of crankshaft and reciprocating components. 


In Figure 2.4, 8; represents the angular position of the crankshaft where the optical 
encoder attaches, 92 through 9, represents the angular position of crankthrows number 


one through three, 9; represents the angular position of the flywheel, and 9. represents the 
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angular position of the load (water brake dynamometer). Each ‘J’ represents the 
equivalent mass polar moment of inertia at each degree of freedom. Equivalent torsional 
stiffness between each degree of freedom is represented by K; through Ks. The material 
damping of the crankshaft between each degree of freedom is represented by Ci2, C23, Ca, 
C4s, and Css. The linear viscous damping effects of each cylinder due to translation of the 
pistons in the cylinders and oil at the journal bearing surfaces are represented by C2, C3, 
and C,. All inertias (J’s) and stiffnesses (K’s) are readily determined by the geometry and 
material properties of the system. Procedures for determining the equivalent inertia and 
stiffness for the model are outlined in Appendix A. The damping for the system, especially 
viscous damping effects, are not so readily determined and will be discussed later when 
comparing measured results to model predicted results. Table 2.1 defines the values of 


inertia (J’s) and stiffness (K’s) for the model [Ref. 8]. 


Table 2.1 Inertia and Stiffness of the model. 


J, = 0.273 + 0.069[1-cos2(0-4n/3)] Ibr-in-sec” 
Ks = 1.304 x 10° Ibpin/rad 











jap Data Required For Model 

Physical dimensions were obtained by measuring the geometry of an equivalent 
crankshaft out of an engine of the same model as shown in Figure 2.3. The material used 
to manufacture the crankshaft is Steel SAE 1548 Modified (Maganese 1.00/1.30 and 


1] 


Carbon 0.43/0.52). | Both the physical dimensions and material properties were checked 
against the manufacturer’s design drawings [Ref. 9] 
Other data required for the model are measured gas torque as a function of angular 


position and load torque 


B. EQUATIONS OF MOTION 


Referring to Figure 2.4, the equations governing the motion of the system are: 


(1) Ji61 + Ci2(6: — 82) + K( 01 — 62) = 0 

(2) (J2+J,., )O2+ C12(O2 — 61) + Ki(82 — 61) + C23(O2 — 3) + K2(02— @3) + C262 = Tiot (t) 
(3) (J3+J;,,. )O3 + C23(O3 — 02) + K2(03 — 2) + C3a(O3 — 04) + K3(03 — 04) + C303 = T2091 (1) 
(4) (Ja+J,.. )Os + C3a(O4— 63) + K3(@4— 03) + Cas(O4 — Os) + Ka(@4 — @s) + C404 = T3o1(t) 
(5) Js0s + Cas(Os— 04) + Ka(Os — 04) + Cs6( Os — Os) + Ks(s — 06) = 0 

(6) J6O6 + Cs6(O6 — Os) + Ks(O6 — Os) = —Thoad (2.2) 


2rec 


The equations of motion in matrix format are: 


7 0 0 0 0 oO1f6, 
Q sett, . 0 0 0 O |/8, 
0 0 J a 0 0 0/8, 
0 0 0 I, +h. 10° 0 || 63 
0 0 0 0 i, OMe. 
0 0 0 0 0 dae, 


C2 —O, 0 0 0 0, 
—C; C2 + Cy a C, —C,, 0 0 0 0, 
0 aCe CRC a Ca 0 0 0, . 
O ¢ sO) Ct Cr ae, C45 0 6, 
O 0 C45 Cyst Cy Ce 0 5 
0 0 0 Cee a 
iG Se 0 0 0 0 0, v 
A, K,+K, ais 0 0 0 0, Dey (t) 
0 aK, K,+K, SK 0 0 O,| | T4.5,(4) 
0 0 -K, K,+K, KK, 0 0, Ds. y7 CO) 
0 0 0 -K, OR ich 0 
0 v 0 0 a i Geal\lle ee ae: 
(2.3) 


The reciprocating component of inertia at each crankthrow, Jorec, Jarec, and Jarec, 
represent the inertia contribution of the piston and connecting rod. The influence of these 
reciprocating components is a maximum when the crankshaft is at the position 
corresponding to mid-stroke (90 degrees after top dead center), and disappears when the 
crank is on top or bottom dead center (0 or 180 degrees after top dead center). The 


relationship of the non-linear reciprocating inertia to crankshaft position is: 


Joree= WR’/2g[ 1-cos(2@)] where: W = Weight of piston and connecting rod (10.5 Ibf) 


Jerec= WR7/2¢f 1-cos2( 0 -27/3)] R = Crankshaft crankthrow eccentricity (2.25 in) 
Varec=WR*/2g[ 1-cos2( 6 -41/3] g = Acceleration of gravity (386 in/sec’) 
8 = Angular position of #1] crankthrow after 
(2.4) TDC (rad) 
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The right hand side of the equations of motion are the applied torques [Ti.yi(t), 
T2-yi(t), and T3,y(t)] into the system and load torque, Tioaa (which is considered constant), 
out of the system. Figure 2.5 shows an example of individual measured cylinder gas 
torque contribution referenced to top dead center of number one cylinder. Since the total 
gas torque into the system has to overcome damping and friction before it exits as useful 
torque (load torque), the gas torques can be multiplied by some factor to account for the 
losses. It is also possible to subtract a constant parasitic torque from each cylinder gas 
torque. Integrating each cylinder torque curve with respect to angular position yields and 
average torque. Dividing the load torque by the sum of the three average cylinder torques 


will yield the mechanical efficiency (nm) to be used. Specifically: 


ow (2.5) 
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Figure 2.5 Individual cylinder gas torque contribution referenced to Nr. | 
cylinder top dead center (TDC). 


c. NUMERICAL SOLUTION OF EQUATIONS OF MOTION 


In order to solve the six second order equations of motion it is first necessary to 


reduce them to a system of twelve first order differential equations. These are: 


| Es) 


0, = A, 
0, = 0, 
0, = 0 
O19 = 0, 
A, = 0, 
A. = 0, 


J,0, — Cy, 9, +C,,0, - K,0, + K,0, = 0 
(J, 20 ao )0, 7 C439, +(C,, + Ci, +C,)@, = C19, 7 K,@, +(K, + K, )6, pa K,@, = L(t) 
(J; ete vas) ep = C349 10 (ORF Bs C3 +C, )6, - C330; a K,6, +(K, 7 K, )0, 7 K,6, = Lo 
(J, cee xo = C459), +(C,, a C34 + €,)A 1, 7 C349, _ K,@, +(K, 1 K,)@, a) K,0; = Loi (t) 
Io _ C5612 + (Cy, 3 C45 )6,, a C4510 - K,9, ar ae K, )6, 7 K,9, = 
EG + C5612 — Coe) + K595 — K5O5 = Tica 

(2.6) 


This initial value problem was solved using a fourth order Runge-Kutta method 


ordinary differential equation solver (ODE45 in MATLAB). Twelve initial conditions are 
required, namely, angular position O(¢=0) and angular velocity 0(¢=0). The system 
repeats after one revolution of the crankshaft (2x rad). A numerical solution for the 


predicted angular position of the first degree of freedom is calculated and presented in 


Section IV. B. 
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HI. EXPERIMENTAL SETUP 


A. THE THREE CYLINDER TWO STROKE DIESEL ENGINE 


il Operating Conditions and Design Specifications 


The test engine used is a Detroit Diesel Series 53 two stroke three cylinder engine. 


The engine design operating characteristics are listed in Table 3.1. 


Table 3.1 Engine Characteristics [Ref. 10] 


Model 
Engine Type 
Number of Cylinders 
Bore and Stroke 
Exhaust Valves per cylinder 
Engine Displacement 
Compression Ratio 
Maximum Power Output 
Full Load Speed 
Peak Torque 
Peak Torque Speed 


Brake Mean Effective Pressure 
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5033-5001N 
In line-2 Cycle-Naturally Aspirated 
3 
3.875 x 4.50 inches 
4 
159 cubic inches 
ZAOal 
92 bhp 
2,800 RPM 
198 ft-lbs 
1,500 RPM 


Ib 
oo 
int 


The manufacturer provides the following estimate of cylinder combustion 
pressures based upon an empirically obtained relationship: 


4,650,000- JHP . _ ; 
= era where: C, = Combustion pressure (psia), Gy 


N =No. of cylinders 
B = Bore (in), S = Stroke (in), 





IHP = indicated horsepower ( 0 


RPM = revolutions per minute of crankshaft. 


= 7.4, = 9./0) 


This formula gives values within + 5% for speeds (1200 RPM to maximum rated 
RPM) and power (continuous to maximum intermittent rated power). [Ref. 11] 


Appendix B [Ref. 12] provides rated power as a function of speed curves and engine 


specification data. 

Z Engine Event Timing 

The major events occurring in the cylinders with respect to crankshaft angle are 
given in Figure 3.1, specifically, intake port, exhaust valve, and fuel injector movements 
with respect to top dead center of the es one cylinder. The firing order of the engine 
(1-3-2) can clearly be seen. The number three cylinder and number one cylinder events 


occur 120 and 240 degrees respectively after the corresponding number one cylinder 


events #3 Cyl intake Ports 


#3 Cyl Exhaust 
#3 Cyl Injector 


#2 Cyt Intake Ports 
#2 Cyl Exhaust 
#2 Cyl Injector 


#1 Cyl Intake Ports 


#1 Cy! Exhaust 
#1 Cyl Injector 


-180 150 i20.~=—--90 60 3C 39 60 90 120 150 180 


Figure 3.1 Engine event timing map referenced to TDC Nr. 1 cylinder 
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B INSTRUMENTATION 

i; Optical Encoder and Mount 

The optical encoder mount holds the casing of the optical encoder while the shaft 
of the encoder is connected to the crankshaft. The mount is in turn attached to the engine 
block. Any torsional displacement between the mount and the engine block will introduce 
errors and hence are to be minimized. The two previous generations of optical encoder 
mounts allowed the encoders to move too much during engine operation. The first 
encoder mount did not provide adequate support for the encoder bearings [Ref. 7] The 
second encoder mount was essentially an aluminum plate mounted on the end of three 
cantilevered beams. When a proximeter was used to measure the displacement of the 
aluminum plate with respect to the engine block, it was determined that the optical 
encoder mount was moving approximately +0.02 degrees with respect to the engine 
block across a wide variation of speeds and torques. This relative displacement, however 
small it seems, corresponds to a 20% error in delta t readings with a 3600 count optical 
encoder, and a 4% error in delta t readings from a 720 count optical encoder. These 
magnitudes of error are unacceptable when changes in delta_t’s between different angular 
positions can vary by as little as three microseconds. A new more rigid optical encoder 
mount was designed and successfully installed. 

Appendix B describes the design of the new optical encoder mount and optical 
encoder installation procedure. Appendix F is the procedure for establishing top dead 
center (TDC) relative to the optical encoder. All events for data collection are referenced 


to TDC of the number one piston. Figure 3.2 shows the new optical encoder mount. 
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Figure 3.2 New optical encoder mount 


2. Pressure Transducers 

These pressure transducers are used to measure gas pressure inside the cylinders 
for cylinder pressure versus crank position data collection. Cylinder pressure is a 
necessary component in the solution of the model predicted response. 

Reference air pressure is a required input for the Engine Cycle Analyzer software. 
These pressure transducers required a reference pressure to properly measure the absolute 
pressure. This reference pressure is the air box pressure or the air pressure the cylinder 
experiences when the intake ports are exposed (see Figure 3.1). A pressure comparator 


connected to the spark plug mounting hole in the side of the engine was used to establish a 
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reference pressure versus engine speed curve. Figure 3.3 is the measured pressure in the 


air box versus engine speed curve. Appendix C describes the procedure used to establish 


the curve. 


Reference Pressure (psi) 


Reference Air Box Pressure vs. Engine Speed 


3.5 


2.9 


500 1000 15:00 2000 2500 
Engine Speed (RPM) 


Figure 3.3 Reference pressure versus speed curve 


3: Modulation Domain Analyzer 


A Hewlett Packard 53310A Modulation Domain Analyzer (MDA) was used to 


measure the time intervals between successive windows of the optical encoder and hence 


determine O(f), O(t), and O(t). This device measures frequency or time interval 


(delta t) between successive events versus time (the modulation domain). The optical 


encoder installed on the engine provides an output signal of a repeating five volt TTL 


2) 


square wave with each rise corresponding to a known angular position (0.1 degrees for a 
3600 count optical encoder and 0.5 degrees for a 720 count optical encoder). With the 
32,000 bit extended memory option, eight revolutions of crankshaft information using the 
3,600 count optical encoder and 44 revolutions of information using the 720 count optical 
encoder (8 revs x 3600 events per rev = 28,800 events and 44 revs x 720 events per rev = 
31,680 events) could be recorded and stored. This information is used to compare 
measured angular position to model predicted angular position in Section IV. C. 

The procedures for measuring delta_t’s from the optical encoder, accessing this 
data from the memory of the MDA, and storing it to facilitate passage to software such as 


MATLAB is discussed in Appendix D. 


4. Data Acquisition Procedure 


Figure 3.4 is a schematic of the test engine and instrumentation. Appendix E 


describes the data acquisition procedure. 


22 


Optical Encoder 


Detriot Diesel (3-53 Series) 3-cyl, 2-stroke 


© @® 
TDC and Crank (Angle Signal TIT 


De eeneasiiiers Amplifiers 


Engine Cycle Analyzer 
Sensor Interface 










— Water brake 
Dynomometer 


Engine Load Signal 






Modulation 
Domain 
Analyzer 


TDC and Crank Angle 


High Speed Data Signals to MDA High Speed Data 
Acquisition Computer 


Acquisition Computer 


Process and store time 
resolved info from 
Optical Encoder 


Process and store 
cylinder pressure info 





Figure 3.4 Test engine and instrumentation 


C. ENGINE CYCLE ANALYZER 

The Engine Cycle Analyzer (ECA) manufactured by SuperFlow 1s a PC based data 
acquisition system which allows pressure versus angular position of the crankshaft data to 
be measured and stored. During data acquisition, the ECA was used to record and store 
pressure data from all three cylinders. The ECA has the ability to monitor | to 999 
revolutions of data from the engine and store the ensemble average of this data to 
memory. When using the 3,600 count optical encoder, 8 revolutions of data were 


recorded and when using the 720 count optical encoder, 44 revolutions of data were 
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recorded. These number of cycles were chosen to match the maximum number of cycles 
that the MDA could store. Revolutions of recorded data corresponded to the data 
collection capability of the MDA. Appendix E discusses operation of the ECA and setup 


files. 
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IV. RESULTS 


A. MEASURED SHAFT SPEED VARIATIONS 

To obtain a sufficient data base, time resolved cylinder pressures and near 
instantaneous shaft speed were acquired for all combination of speeds and torques shown 
in Table 4.1. These represent low, medium, and high torque and low, medium, and high 
speed. For each condition (e. g. 1,000 RPM, 135 FT-LBS) data were acquired with both 


the 720 count and the 3,600 count optical encoder. 


Table 4.1 Data mins 


SPEED MEASURED LOAD TORQUE 
RPM FT-LBS) 


| ooo | 22) | 80 | 385 









Only the result of the 1,000 rpm and 135 ft-lb mun utilizing the 720 count optical 
encoder was chosen for analysis. Raw data from the Modulation Domain Analyzer was 
phase locked and ensemble averaged to represent one typical cycle (0 to 360 degrees after 
TDC). Appendix G is the MATLAB code which phase locks and ensemble averages the 
raw data. Figure 4.1 is a plot of the phase locked, ensemble average data from the MDA 


for 44 cycles taken at 1,000 RPM and 135 FT-LBS torque. 
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Figure 4.1 MDA phase locked, preemble average of raw data taken at 1,000 rpm 

and 135 ft-lbs. 

Figure 4.1 shows on the ordinate the time between successive windows (delta _t) 
for the optical encoder versus time. This corresponds to 360 degrees of data. If the 
crankshaft were turning with a constant angular velocity, the data shown in Figure 4.1 
would be a horizontal line. Any value below the average delta t indicates a speed-up of 
the crankshaft (it takes a shorter time for the next window to arrive). Any value above the 
average delta t indicates a slow down of the crankshaft (it takes more time for the next 
window to arrive). The average angular velocity for the 1,000 rpm, 135 ft-lb run is 
107.268 rad/sec. Appendix G is the MATLAB code that computes average angular 


velocity. The average time between pulses 1s 81.35 micro-séconds. 
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These data can be presented in any of several ways. The reciprocal of delta t 
yields the shaft speed. However, these data tend to be quite noisy. Therefore, the data 
will be processed and presented using the deviation in angle from a constant speed shaft. 

With the ensemble averaged delta_t’s and the average angular velocity computed, 
the measured response at the free end of the crankshaft can be calculated and plotted. The 


equation to calculate is: 


O eased = pom ~ eycte * tmearred Where: Oy, = angular position of optical encoder 
window relative to TDC 
® ce — average angular velocity of crankshaft 


(4.1) t easured = time to arrival of optical encoder 
window 


Plotting the angle of the crankshaft as opposed to the angular velocity response reduces 
noise. Taking the integral of the angular velocity produces results that are less sensitive to 
noise. If the crankshaft was operating at a constant angular velocity, the response plot 
would be a horizontal line. Deviation from constant velocity indicates that the crankshaft 
angular position is leading or lagging a constant velocity phasor. Appendix G is the 
MATLAB code that computes measured crankshaft position deviation at the free end of 
the crankshaft. Figure 4.2 is a plot of measured crankshaft position deviation from a 


constant speed shaft. 
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Figure 4.2 Measured Response of crankshaft at free end 
B. 


TORSIONAL VIBRATION MODEL PREDICTED RESULTS 
Utilizing the measured pressures and calculated inertia/stiffmess from the model, 
the equations of motion were solved numerically utilizing a fourth order Runge-Kutta 
method. The initial conditions were: H(t = 0) = 0 rad and A(t = 0) = 108.9 rad/sec for all 
degrees-of-freedom. The material damping was chosen as 0.002 and viscous damping 
0.02. The non-linear reciprocating inertia terms were updated every 0.5 degrees or 720 
times per revolution. The system of 12 first order differential equations was integrated 
from tinitiat = 0 tO tinal = 0.058 seconds. There were 720 linearly spaced time steps between 
tinitiat ANd tfna. These time steps were chosen to match the ECA data updates. At the 
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completion of each time step, the response of 8; was recorded. Appendix G is the 


MATLAB code that computes predicted response. Figure 4.3 is a plot of predicted 


crankshaft position deviation from a constant speed shaft. 
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Figure 4.3. Model predicted response of crankshaft at free end 


COMPARISON OF PREDICTED TO MEASURED RESPONSE 


Figure 4.4 shows a comparison of predicted and measured response versus time. 


Several similarities are obvious. First, the slope reversals occur at nearly the same phase 


and time. Secondly, the slopes are generally the same. Measured and predicted responses 


show good qualitative agreement. 
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Figure 4.4 Measured and predicted responses of crankshaft at free end 


To better understand these curves it is necessary to consider the driving torque. 
Figure 2.5 shows that the maximum gas torque for the number one cylinder occurs 20 
degrees after TDC. With the crankshaft rotating at 107.3 rad/sec (1024 rpm), 20 degrees 
after TDC of number one cylinder occur approximately 0.0033 seconds into the cycle. 
The system speeds up (positive slope) until the applied torque is zero at approximately 80 
degrees (0.013 seconds) after TDC. When the applied torque is zero the slope of the 
predicted response reverses. This pattern continues for the firings of cylinders three and 
two. All three distinct increases in angular velocity of the crankshaft associated with 


cylinders firing are clearly present. 


The measured response with all of the stated similarities still has one large 
difference. At the beginning of the cycle, the slope decreases and continues this pattern 
until 0.025 seconds after TDC. After 0.025 seconds, the angular velocity increases until 
the cycle repeats. Examining the measured response curve indicates that the applied 
torque into the system is lower at the beginning of the cycle and then increases toward the 
end of the cycle allowing repeatability. Figure 2.5, the individual cylinder gas torque 
curves, show low gas torque for the number one cylinder, slightly increased gas torque for 
number three cylinder, and finally the greatest gas torque contribution by number two 
cylinder. This would appear to explain the shape of the measured response curve, low 
pressure in number one cylinder. 

To validate the measurement, pressure transducers were interchanged between 


cylinders. The same pressure curves were obtained. 





V. CONCLUSIONS AND RECOMMENDATIONS 


A. SUMMARY 

A new rigid optical encoder mount was designed and installed on the engine. A 
torsional vibration model was developed which included the crankshaft, major 
reciprocating components, and external loads. Equations of motion describing the 
dynamics of the system were derived. The Modulation Domain Analyzer and Engine 
Cycle Analyzer were used to gather time resolved information from the engine. Baseline 
data were obtained at eighteen different speed and torque combinations utilizing a 720 and 
3,600 count optical encoder. The predicted response and the measured response were 


compared. 


Cc. CONCLUSIONS 

The predicted and measured response at the free end of the crankshaft show good 
qualitative agreement. The predicted response shows the expected increase and decrease 
in instantaneous angular velocity as a result of cylinders firing. The torsional vibration 
model with non-linear inertia and linear stiffness and damping appears to describe the 
dynamics of the diesel engine well. Specifically, the six degrees-of-freedom is sufficient to 
capture the measured torsional dynamics. Measured response appears to be different due 


to actual low cylinder pressure in number one cylinder. 
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C. RECOMMENDATIONS 

The applied torque to the shaft should be calculated using a parasitic torque, where 
a constant magnitude is subtracted from each cylinder instead of multiplying each gas 
torque by the same mechanical efficiency. This should greatly improve the agreement 
between the model and the measurements. 

Other contributions of inertia from the camshaft, counter-balance shafting, gears, 
and exhaust valves should be explored for significance. 

More model validation can be obtained by instrumenting the flywheel with a 
proximeter. This adds instantaneous angular velocity information from the fifth degree of 
freedom. The torsional twisting between the first and fifth degrees of freedom may also ~ 
hold information of interest. 

The low pressure in the number one cylinder should be investigated and corrected. 
Once corrected, another 1,000 rpm, 135 ft-lb torque data run should be made. The 
measured response curve can be compared to the measured response from this research. 
It is expected that the results will be similar to the model predicted response. 

The next step should be pressure prediction based upon instantaneous shaft speed 
measurement. The measured response curve could also be used to localize cylinder and 


engine faults. 
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APPENDIX A 


INERTIA AND STIFFNESS FOR THE TORSIONAL ENGINE MODEL 


To calculate the equivalent inertia and stiffness of the model, Reference 1 was 
utilized. In general, the complex system must be converted into a simpler system of 
equivalent masses and equivalent shafts. The manufacturer’s design drawings [Ref. 9] of 
the crankshaft were utilized to obtain the geometry and material properties. 

MASS POLAR MOMENT OF INERTIA OF CRANKSHAFT AND RECIPROCATING 
COMPONENTS: 

Il. Each equivalent mass polar moment of inertia of the crankthrows include all inertia 
contributions from the mid-point of a journal bearing to their left and to the their nght, 
Figures 2.3 and 2.4 germane. The inertia from the mid-point of the front journal bearing 
surface to the end of the crankshaft is J;. It includes the inertia contribution of the optical 
encoder crankshaft adaptor. The inertia from the mid-point of the rear journal bearing 
surface to the back of the crankshaft including the flywheel is J;. The equivalent inertia of 
the water brake dynamometer is Jc. The value of Js was obtained from the manufacture of 
the water brake [Ref. 13]. The reciprocating inertia from the pistons and connecting rods 
is discussed later. 

D, In general, to compute the equivalent inertia of the system, the crankshaft must be 
broken down into simpler components. The following are the component inertias which 


must be summed to obtain the equivalent inertia at each degree of freedom. 


a. Crankshaft journals (solid): 
W ==*D?* A* p (Ibf) 
= r A p 
2 
Kos a (in’) 


* 2 
J= dees (Ibf-in-sec’) 
& 





where: W = weight of journal in Ibf 
K = radius of gyration of journal in inches 
D, = diameter of journal in inches 
A = length of journal in inches 
p = specific weight of material = 0.283 Ibf/in° 
g = 386 in/sec’ 


b. Crankpins (solid): 
2 * 2 
ee aa epeiecec’) 
1.735 | 8 
where: D2 = diameter of crankpin in inches 
B = length of crankpin in inches 
R = distance of center of gravity of crankpin 
from the crankshaft axis of rotation, in inches 





c. Crankwebs: 
W * K’ 


= 
a 


where: W = total weight of component in lbf 
K = radius of gyration in inches 
g = 386 in/sec’ 
note: Each crankweb must be broken down into small geometric 
shapes. The parallel axis theorem is utilized to transfer 
the equivalent inertia to the crankshaft axis of rotation. 
Table 10.13 [Ref. 1] contains radi of gyration (K) and 
weights (W) for 28 geometric shapes. 
d. Reciprocating inertia: 
2 
Po ae) (Ibf-in-sec’) 
2 * g[1 — cos(28)] 
where: W = weight of piston and connecting rod = 10.5 lbf 
R = distance of center of gravity of crankpin 
from the crankshaft axis of rotation, in inches 
g = 386 in/sec” 
8 = angular position of crankthrow after TDC 





(Ibf-in-sec’) 
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EQUIVALENT TORSIONAL STIFFNESS OF CRANKSHAFT: 

1. Each equivalent torsional stiffness was calculated from one degree of freedom to 
the next. Torsional stiffness between crankthrows is from the mid-point of the crankpin of 
one cylinder to the mid-point of the crankpin of the next cylinder. 

Z. In general, to compute the equivalent torsional stiffness of the system, the 
crankshaft must be broken down into equivalent lengths of shafting. Component 
equivalent lengths must be summed to get equivalent lengths between each degree of 


freedom. 


a. For journals (solid), crankpins (solid), and crankwebs: 


* * os 
L.= p* GAN D;\ cant IDs 4 ead (in) 
D; Ds T*W 


one journal one pin two crankwebs 


where: D = diameter of equivalent shaft in inches 
A = length of journal in inches 
D, = diameter of crank journal bearing in inches 
D, = diameter of crankpin 
R = distance of center of gravity of crankpin 
from the crankshaft axis of rotation, in inches 
T = thickness of crankweb 1n inches 
W = width of crankweb in inches 
4 
Kee Bin/ad) 
Bo ee 
where: G = modulus of rigidity = 12.0 X 10° (Ibf/in’ ) 


b. For solid shafts: 
* [)* * 
Re ela i iced) 
Bye 
c. To calculate the equivalent stiffness from crankpin to crankpin, calculate the 
stiffness as shown above (Kmain). Main indicates main journal bearing surface to 
main journal bearing surface. Koin-to-pin = (2*Kmain) + (2*Kmain) 
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d. Figure 11.5 [Ref. 1] contains some typical crankshaft elements and equivalent 
lengths. 
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APPENDIX B 
OPTICAL ENCODER AND OPTICAL ENCODER MOUNT INSTALLATION 
PROCEDURE 

In order to obtain accurate readings from the optical encoder, the following 
procedure should be followed when installing the mount and or changing optical encoders: 
Installing mount: 
i Place the crankshaft to optical encoder stainless steel adapter on the end of the 
crankshaft and firmly seat. Rotate the adapter until the two encoder shaft set screws are 
pointing straight up and down or vertical to the deck. Tighten the three set screws 
securely against the crankshaft with an allen wrench. 
5, Place the smaller of the two aluminum optical encoder mounts against the engine 
block and then place the larger one on top oriented such that the encoder will mate. 
Obtain the optical encoder mount bolts and hand tighten the bolts until the mounts are 
flush with the block. Figure 3.2 shows this orientation. 
3. Obtain the optical encoder alignment jig. Place it into the countersunk encoder 
receptor. The purpose of the jig is to align the mount so the crankshaft rotates 
concentrically with the optical encoder shaft. Carefully tighten down all three optical 
encoder mounting bolts. Ensure that the jig 1s able to move in and out of the encoder 
adapter hole and the mounting lip without binding. It may be necessary to shim the 
encoder mount. If all three bolts are securely fastened and the jig moves in and out freely 


from the receptor, the optical encoder mount is aligned. 


oy 


Installing optical encoder: 

l. Place optical encoder into the receptor and install the three small bolts which retain 
it. Tighten the bolts with an allen wrench. 

pA On the bottom of the optical encoder there is a window which allows access to the 
encoder shaft set screws. Tighten the one set screw securely against the optical encoder 
shaft with an allen wrench. It may be necessary to utilize a flashlight and mirror to 
accomplish this. 

3. Remove the access plate over the flywheel. Rotate the engine crankshaft 
clockwise with a pry bar 180 degrees. This allows the other the set screw to be accessed. 
4. Tighten the remaining set screw securely against the optical encoder shaft with an 


allen wrench. 


The optical encoder mount and optical encoder are now aligned. Ensure that TDC 
relative to the optical encoder is established prior to engine operation, Appendix G 


germane. 
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APPENDIX C 


ESTABLISHING REFERENCE PRESSURE 


In order to obtain the most accurate cylinder pressure readings, the ECA requires a 
reference pressure input. This appendix describes the procedure for establishing reference 
pressure versus speed. It is obtained at normal operating temperature of the engine. 

An MKS Instruments Inc. Type 698 differential high-accuracy pressure transducer 
in conjunction with the MKS Instruments Inc. Signal Contitioner Type 270B was used to 
measure engine air box pressure relative to atmospheric pressure. The range of the 
differential pressure transducer was 1,000 torr (20 psi), sufficient to measure the air box 
pressure. The Air box pressure was obtained by removing a cold weather spark plug on 
the side of the engine and using the threaded hole for a air box pressure tap. The other 
side of the differential pressure transducer was open to the atmosphere. The engine was 
operated from 700 to 2,500 rpm in increments of 200 rpm. At every increment of speed, 
the air box pressure was recorded. A plot of reference air box pressure versus engine 
speed was constructed, Figure 3.3 germane. When measuring cylinder firing pressures, 
the reference air box pressure for the speed of engine operation should be input into the 
ECA set up file. 

The following picture shows the set up of the pressure transducer, signal 


conditioner, and associated piping to the engine. 


4] 





Figure C.1 Instrumentation utilized for establishing reference air box pressure 


versus engine speed. 
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APPENDIX D 


MODULATION DOMAIN ANALYZER OPERATING PROCEDURE 


This appendix describes how to configure and operate the modulation domain 


analyzer (MDA) to record time resolved information from the optical encoder during 


engine operation. The procedure is: 


le 


a 


eo 


ale 
ze 


D:: 
14. 
15. 


16. 


17. 
M7. 
18. 
Ze 


Ensure both the MDA and data acquisition computer are turned off. 

Connect one end of the IEEE-488 cable (National Instruments 763061-02 Rev C, 
Type-X2, Length 2.1 meter) to the HPIB on the back side of the MDA. Connect 
the other end of this cable to the National Instruments GPIB board interface on 
back of the data acquisition computer. 

Connect the TDC output on front of the ECA sensor interface to the EXT Arm 
input on front of the MDA with a BNC connector. 

Connect the CA output on front of the ECA sensor interface to the Channel A 
input on front of the MDA with a BNC connector. 

Turn the power to the MDA on. 

Depress the ‘function and input’ menu key. 

Select the ‘time interval’ soft key. 

Select ‘common’ A and B line soft key. 

Select ‘rising edge A and Rising edge B’ soft key. 

Depress the ‘trigger’ menu key. 

Select ‘triggered’ soft key. 

Select ‘mean’ and ‘1/mean’ on the entry panel to display the mean delta_t and 
mean frequency. 

Start the engine in accordance with Appendix A [Ref. 1]. 

Depress the ‘auto scale’ menu key. 

The MDA is now displaying delta_t on the y-axis and time on the x-axis. Continue 
with step 16 to make the MDA ready for remote operation and data acquisition. 
Turn on the power to the data acquisition computer and go to the HP Applications 
icon and enter. | 

For remote operation, depress the utility’ menu key on the MDA. 

Select ‘addressed’ under the HP-IB/Print soft key. 

Select address number 5 (five). 

The MDA is now in remote operation mode. In order to print the screen or 
operate the MDA locally, the ‘TALK ONLY’ mode must be enabled. 
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The following section discusses how to access the memory of the MDA, record 


the measured delta _t’s, and store this information with the data acquisition computer. 


el 


le 


With the MDA properly configured, select the HP Applications icon and go to the 
HP53305A Phase Analyzer software. Connect ‘CA’ lead to ‘B’ input on MDA. 
Under the setup menu, select instrument setup and configure the following: 

a. Input configure: B, + slope, DC coupling, and 50 ohm impedance. 

b. Arming: EXT, TTL, and + slope. 

c. Interface: National GPIB, 0 board index, and address 5. 

Under the setup menu, select measure/view setup and configure the following: 

a. Parameter: Frequency and Phase. 

b. Number of samples: Enter the number of sample points desired for the run but 
do not exceed 32,000. Remember, the 720 count optical encoder is 720 data 
samples/revolution and the 3,600 optical encoder is 3,600 samples/revolution. 

c. Measurement time: Enter the approximate expected delta_t’s expected for the 
speed you are taking data at. (1.e. at 1000rpm, with the 720 optical encoder, 

min , rev ,, 00 sec 

1000rev 720windows min — 





enter 83.3 microseconds or 


. Number of averages: one. 

. Pacing: auto. 

. Input: B. 

Clock/carrier frequency: calculated. 

. No filter. 

Under the setup menu, select phase units then select ‘unit interval’. 

The MDA is now ready to collect TDC triggered data from the optical encoder. 
When ready to start collecting data, go to the measure menu and select ‘phase 
deviation’. As soon as you release the mouse button, the MDA records data 
Once the data is recorded, the MDA will ask you what to name the plot. This 
is unimportant, enter anything. 

Save the data to the hard drive or zip drive as a M-File. 

With MATLAB, remove the header and footer information which leaves two 
columns of numbers. The length of the array is the number of sample points you 
previously selected during setup of the software. 

The data file is now ready for use by MATLAB code, Appendix H germane. 


ar (2ty te (ee 
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Figure D.1 MDA and high speed data acquisition computer. An oscilloscope sets 
on top of the MDA. 





APPENDIX E 


DATA ACQUISITION PROCEDURE 


This appendix describes the data acquisition procedure for comparing model 


predicted response to measured response at the free end of the crankshaft utilizing the 720 


count optical encoder. The procedure is: 
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Align the engine for operation and turn the charge amplifiers on. Ensure the 
proper scale and sensitivity is set on charge amplifiers. This information is on the 
Kistler pressure transducer calibration sheets. 

Setup the ECA setup file used for data acquisition in accordance with Appendix F. 
Align the MDA and data acquisition computer in accordance with Appendix D. 
Increase the engine and load to the desired set point. 

On the ECA, select data taking and statistics under the main menu. 

Select 720DD353.SET under the data taking and statistics menu. 720DD353.SET 
is the preset for monitoring and recording all three cylinder pressures utilizing the 
720 count optical encoder. DD353PRE.SET is the setup file for monitoring and 
recording all three cylinder pressures utilizing the 3,600 count optical encoder. 
Enter the correct reference pressure from Figure 3.3 into the J.1 block and enter 
the number of cycles to be samples (should correspond with the MDA). 

Place the charge amplifiers in the ‘operate’ mode. 

Enter ‘G’ for graphic monitoring and ensure three pressure versus theta curves 
appear. 

The ECA is now ready to record data. 

The ECA will temporarily turn on and off the output crank angle mm TDC signals 
while it is acquiring data. Therefore, the MDA data must be taken first. 

Record the speed and torque of the engine by depressing ‘Record Data’ on the 
SuperFlow Dynamometer Control Console. 

Take the MDA measured data in accordance with Appendix D. 

Immediately after the MDA data acquisition data is taken, depress “D” for data 
taking on the ECA keyboard. Follow on screen instructions. 

Store each individual cylinder pressure to the hard drive. 

Data acquisition is complete. Steps 12 through 15 above require practice and 
should be accomplished as fast as possible to avoid any deviation in speed and 
torque from the engine. 
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Figure E.2 SuperFlow Dynamometer Control Console 
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APPENDIX F 
ESTABLISHING TOP DEAD CENTER RELATIVE TO THE OPTICAL 
ENCODER 
In order to get accurate readings from the optical encoder, the following 
procedure should be followed to correctly position Top Dead Center (TDC) of the 
number one cylinder relative to the optical encoder. TDC is the reference against which 
all angular measurements are made. It is the point in the cycle of a cylinder when the 
volume is the least. The number one cylinder is the cylinder for this reference. The 


procedure outlined in Appendix G of Reference 1 is unchanged with the following 


exceptions. 

1. To rotate the engine, remove the flywheel access cover on the back of the engine. 

a A pry bar, as the one depicted, can be utilized to rotate the engine. 

3 Once mechanical TDC of the number one cylinder is established, loosen the three 
retaining bolts that secure the optical encoder to the mount. 

4. Rotate the optical encoder until the TDC signal on the ECA sensor interface 
illuminates. 

=) Carefully tighten the retaining bolts on the optical encoder base. 
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Figure F.1 Flywheel access cover removed with pry bar in place. 


APPENDIX G 


MATLAB CODE USED FOR CYCLE ANALYSIS 


This Appendix contains all the different MATLAB M-files used for this research. 
Data from both the ECA and MDA require that the footer and header information be 
removed prior to loading into an M-file. The ECA records data pressure from -180 
degrees to 180 degrees relative to TDC. Before loading the ECA data into an M-file for 
analysis, it must first be rearranged and saved from 0 degrees to 360 degrees after TDC. 
This procedure must be followed for each cylinder pressure data. 

The following code can be utilized to: 

a. Plot pressure versus crank angle for all three cylinders. 

b. Plot gas torque versus crank angle for all three cylinders. 

c. Compute average gas torque contribution from each cylinder 

d. Compute the mechanical efficiency based upon a torque balance. 

e. Compute and plot the phase lock, ensemble average of MDA raw data. 

f. Compute the average angular velocity through each revolution from MDA data. 

g. Solve the equations of motion of the system and plot a predicted response at 
the free end of the crankshaft. 

h. Compute and plot the measured response at the free end of the crankshaft 
utilizing MDA phase lock, ensemble average data. 

[. Define the equations of motion as a system of 12 first order differential 


equations. 
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% This code computes and plots the Individual cylinder pressures 

“% Referenced to Nr. 1 cylinder TDC 

% Filename: cylpress.m 

load d15kr135.m % ECA cylinder #1 pressure data 

load e15kr135.m % ECA cylinder #2 pressure data 

load f15kr135.m % ECA cylinder #3 pressure data 

% 

“PLOTS PRESSURE OF EACH INDIVIDUAL CYLINDER VS DEGREES AFTER TDC OF NRI1 CYL 
% 

theta=1:1:720; 

plot(theta/2,(d15kr135(theta)* 14.5),'x' theta/2,(e15kr135(theta)* 14.5),'o',theta/2,(f15kr135(theta)* 14.5),'- 
<1) 

title(‘Individual Cylinder Pressure referenced to TDC of Nr. 1 cylinder’) 

xlabel(‘Degrees after TDC of Nr. 1 Cylinder’) 

ylabel(‘Pressure (PSI)') 

legend(‘Nr 1 cyl pressure’,"Nr 2 cyl pressure'’,"Nr 3 cyl pressure’) 

axis([-inf 361 0 1250]) 


% This code computes the individual torque contribution of each 

% individual cylinder referenced to TDC of Nr.1 Cylinder and the 

% recipricating inertia. It can also be manipulated to plot 

% pressure vs. theta curves 

% Filename: torque.m 

% 

load dlkr135.m % ECA cylinder #1 pressure data 

load elkr135.m % ECA cylinder #2 pressure data 

load f1kr135.m % ECA cylinder #3 pressure data 

r=2.25;W=10.5;g=386; “oconstants 

for theta=1:720 
Tlcyl(theta)=(d1kr135(theta)* 14.5*r*sin(theta*.008726646)*(pi*r’2))/12; % FT-LB 
T2cyl(theta)=(e 1kr135(theta)*14.5*r*sin((theta*.008726646)-(4*pi/3))*(pi*r*2))/12;% FT-LB 
T3cyl(theta)=(flkr135(theta)* 14.5*r*sin((theta* .008726646)-(2*pi/3))*(pi*r’2))/12:% FT-LB 
j2rec(theta)=(W*r°2/(2* g))*(1-cos(2*theta*.008726646)); %lb*in*sec*2 
j3rec(theta)=(W*r42/(2* g))*(1-cos(2*(theta*.008726646)-(2*pi/3))); Ylb*in*sec*2 
j4rec(theta)=(W*r*2/(2*g))*(1-cos(2*(theta*.008726646)-(4*pi/3))); Ylb*in*sec*2 


% plot(theta/2,(T lcyl+T2cyl+T3cyl),’kx’) 

% title('Net Torque Input Referenced to TDC of Nr. 1 Cylinder’) 
% Xlabel(‘Degrees after TDC of Nr. 1 Cylinder’) 

% ylabel('Gas Torque (FT-LB)') 

% hold on 

end 

figure(1) 

theta=linspace(1,720,720); 

plot(theta/2,T Icyl,'X',theta/2,T2cyl,'O',theta/2,T3cyl,".') 
title(Individual Cylinder Torque input Referenced to TDC of Nr. 1 Cylinder’) 
xlabel('Gas Torque (FT-LB)’) 

ylabel('Degrees after TDC of NR. 1 Cylinder’) 


a2 


legend(‘Nr 1 Cyl Torque','Nr 2 Cyl Torque',"Nr 3 Cyl Torque’) 

figure(2) 

plot(theta/2.j2rec.'X' theta/2 j3rec,'O' theta/2,j4rec,'.') 

title(‘Individual Reciprocating Crankthrow Inertia Referenced to TDC of Nr.1 Cylinder') 

xlabel(‘Inertia (lb*in*sec’%2') 

ylabel('Degrees after TDC of Nr. 1 Cylinder’) 

legend(‘Nr 1 Cyl Reciprocating Inertia’,"Nr 2 Cyl Reciprocating Inertia’,"Nr 3 Cyl Reciprocating Inertia’) 


% 

% THIS SECTION COMPUTES THE AVERAGE TORQUE CONTRIBUTION OF EACH CYLINDER 
% 

tic 

x=linspace(0,(2* pi), 720): 

Tlcyl_avg=trapz(x.T lcyl)/(2* pi) 

T2cyl_avg=trapz(x,T2cyl)/(2*pi) 

T3cyl_avg=trapz(x,T3cyl)/(2*pi) 

Avg Torque_input=(Tlcyl_avg+T2cyl_avg+T3cyl_avg)/3 

Torque_out_to Torque_in=134/(TIcyl_avg+T2cyl_avg+T3cyl_avg) 


% Jim Hudson Thesis file 
% 720 OE, SAMPLE RUN 1000 RPM, 135 FT-LB TORQUE 
% 16 NOV 1997 Filename: ensemble.m 
% 31680 sample points from MDA or 44 revs of INFO 
% THIS CODE PHASE LOCKS AND ENSEMBLE AVERAGES 
% THE OUTPUT REPRESENTS ONE REVOLUTION OF INFO 
load 71kr135.m % This is where the 31680 times are located 
t=X71kr135(:,1); “% Assigns vector "t" to each time sample 
for b=1:3 1680 

tt(b)=t(b+1)-t(b);  % Computes the del_T's and assigns 
end % them to Vector "tt" 
tt=(reshape(tt,720,44))'; % Phase locks and adds 

=sum(tt)/44; % Averages the phase lock 

position=linspace(0,(2*pi), 720); 
timer=0; % Initialize a counter 
for e=1:720 

time=[e;t(e+1);tttt(e)]; 

plot(timer,tttt(e),"kX") %PLOTS TIME VS. DEL_T 

timer=timer+tttt(e); 

hold on 
end 
gtext('l000 RPM & 135 Ft-Lbs, 720 Optical Encoder’) 
title(MDA DATA PHASE LOCK ENSEMBLE AVERAGED") 
ylabel(‘Delta_t btwn windows (Seconds)') 
xlabel("Time of sample, after TDC (SECONDS)') 


“% This code computes the average angular velocity through each 
% of the 44 revolutions of of MDA data. Additionally, it gives 
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% the average angular velocity the MDA saw. 
% filename: omegacyc.m 
% 
load 71kr135.m % This is where the MDA raw data is stored 
for a=1:31681 
t(a)=X7 1kr135(a); 


end 

for n=1:44 
x(n)=720*n: 

end 


omega _ zero=(2*pi)/(t(720)-t(1)); 

plot(0O.omega_zero,'kX") 

hold on 

sum=omega_ Zero; 

for c=1:43 
omega_avg(c)=(2*pi)/(t(x(c+1))-t(x(c)+1)); 
sum=sumtomega avg(c); 

end 

c=linspace(1,43,43); 

average=sum/44 % This is the average angular velocity 

plot(c,omega_avg,'kX') 

titleAverage Angular Velocity Through Each of the 44 Revolutions’) 

xlabel('Cycle Number') 

ylabel(‘Average Angular Velocity (rad/sec)') 


% THIS CODE PLOTS CONCURRENTLY THE MEASURED AND PREDICTED RESPONSES 
% for a 720 count optical encoder. The measured and predicted 
% response can be plotted independently by simply commenting out 
% one orthe other. Filename: measpred.m 
% 31680 sample points from MDA are located in file: 71KR135.M 
% ECA cylinder pressure information is in files: 
% dikr135.m cylinder #1 pressure data 
% €1kr135.m cylinder #2 pressure data 
“% f1kr135.m cylinder #3 pressure data 
% 
% This first section plots the measured response 
% 
load 71kr135.m % THIS IS WHERE THE 31680 TIMES ARE LOCATED 
t=X71kr135(:,1); % ASSIGNS VECTOR "t" TO EACH TIME SAMPLE 
for b=1:31680 

tt(b)=t(b+1)-t(b); % COMPUTES THE DEL_T'S AND ASSIGNS THEM TO 
end % VECTOR "tt" 
tt=(reshape(tt,720,44))':% This phase locks one cycle 
tttt=sum(tt)/44; % This ensemble averages the phases 
position=linspace(0,(2*pi),720); % The known positions of the O.E. windows 
timer=0; 
for e=1:720 

plot(timer,(position(e)-(107.278*timer)),'rO') %107.278 is omega_bar 

timer=timer+tttt(e); 

hold on 
end 
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gtext(‘Measured Response’) 
hold on 
% 
% This section is a SIX SECOND ORDER SUMULTANEOUS EQUATION ODE SOLVER 
% This section calls "eqns.m" which defines the system of 2nd order ode's 
% 
global Tload plcyl p2cyl p3cyl T1cyl T2cyl T3cyl ... 
j2rec j3rec j4rec jl j2 j3 j4j5 j6 k1 2 k3 k4... 
k5 cl2 c23 c34 c45 c56 c2 c3 c4 “share the variables with other m-files 
ic=[(0 00000 108.909 108.909... 
108.909 108.909 108.909 108.909]; %edefine IC's [x(t=0) & xdot(t=o)] 
load d1kr135.m % cylinder #1 pressure data from the ECA 
load elkr135.m % cylinder #2 pressure data from the ECA 
load flkr135.m % cylinder #3 pressure data from the ECA 
Tload=1608; % lb-in (134 FT-LB as measured by the dynamometer) 
z=linspace(0,.057692066,721); % divide .058 secs into 720 evenly 
% spaced time increments 

r=2.25; % Crankshaft eccentricity (in) 
W=10.5; % The weight of the reciprocating components 

% (piston and connecting rod) (Ibf) 
g=386; % The acceleration of gravity (in/sec’2) 
mech_ eff=.6175; % This is the scale factor for the Gas Torque in 

% mech_eff=Tload/(T lavg+T2avg+T3avg) 

“% M-file Torque.m computes the gas torques 
omega_bar=109.9; %This is the mean speed of the crankshaft used 

“by the plot command for "Response Deviation" 
jl=.0291; %lb*in*sec*2/rad 
j2=.2730, %“%lb*in*sec’2/rad 
j3=.1710; “%lb*in*sec’*2/rad 
j4=.2730; %“%lb*in*sec’2/rad 
j5=6.194; %lb*in*sec*2/rad 
j6=.2800; %lb*in*sec’2/rad 
k1=3.470e6; %“olb*in/rad 
k2=9.551le6; “%olb*in/rad 
k3=9.551le6; “lb*in/rad 
k4=13.495e6; “%olb*in/rad 
k5=1.304e6; “lb*in/rad 
c12=.002; %lb*in*sec/rad 
c23=.002; %lb*in*sec/rad 
c34=.002; %lb*in*sec/rad 
c45=.002; %lb*in*sec/rad 
c56=.002; “%lb*in*sec/rad 
c2=.02; %lb*in*sec/rad 
c3=.02: %lb*in*sec/rad 
c4=.02; “%lb*in*sec/rad 
X=I1C; Yolnitialize the "x" vector to compute Tcyl's and Jrec's 
for theta=1:720 
Tlcyl=mech_ eff*d1kr135(theta)* 14.5*r*sin(theta*.008726646)*(pi*r%2);“olb-in 
T2cyl=mech_eff*e1kr135(theta)* 14.5*r*sin((theta*.008726646)-(4*pi/3))*(pi*r*2):%lb-in 
T3cyl=mech_eff*flkr135(theta)* 14.5*r*sin((theta*.008726646)-(2 *pi/3))*(pi*r*2):“elb-in 
j2rec=(W*1r*2/(2* g))*(1-cos(2 *theta*.008726646)); “olb*in* sec’*2 
j3rec=(W*1r*2/(2* g))*(1-cos(2 *(theta*.008726646)-(2*pi/3)));  Yelb*in*sec%2 
j4rec=(W*1r*2/(2* g))*(1-cos(2 *(theta*.008726646)-(4*pi/3))); Yolb*in*sec%2 
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[T.x]=ode45(‘eqns',z(theta).z(thetat 1),ic); 
ic=x(length(T),:); % REINITALIZE THE IC'S FROM THE PREVIOUS ITERATION 
plot(T(length(T)),x(length(T), 1)-(omega_bar*T(length(T))).’bX’) 
end 
title(‘Crank Position (PREDICTED/MEASURED) Deviation from Mean‘) 
xlabel('Time (sec)') 
ylabel(‘Angular Position (rad)') 
axis({O .0573 -.027 .01]) “% Change this for the appropriate desired view 
gtext(‘Predicted Response’) 


% This code defines the 12 first order equations to be differentiated. 

% FILENAME: EQNS.M 

function xdot=firsteqn(t,x,Tload) 

global Tload Tlcyl T2cyl T3cyl j2rec j3rec j4rec jl j2 j3... 
j4 j5 J6 k1 k2 k3 k4.kS cl2 c23 c34 c45 c56 c2 c3 c4 

% 

“% These are 12 first order equations which define the 

% Original 6 second order equations of motion!!!! 

% 

xdo(1)=x(7); 

xdo(2)=x(8); 

xdo(3)=x(9): 

xdo(4)=x(10); 

xdo(5)=x(11); 

xdo(6)=x(12); 

xdo(7)=((c12*x(8))-(c12*x(7))+(k1*x(2))-(k1*x(1)))41; 

xdo(8)=(T lcyl+(c23*x(9))-((c23+¢12+¢2)*x(8))+(c12*x(7))+(k2*x(3))... 
-((k1+k2)*x(2))+(k1*x(1)))/G2+j2rec); 

xdo(9)=(T2cyl+(c34*x(10))-((c34+¢23+¢3)*x(9))+(C23 *x(8))+(k3 *x(4))... 
-((k2+k3)*x(3))+(k2*x(2)))/Q3+j3 rec); 

xdo(10)=(T3cyl+(c45*x(11))-((c45+¢34+c4)*x(10))+(c34*x(9))+(k4*x(5))... 
~((K3+k4)*x(4))+(k3*x(3)))/G4+)4rec); 

xdo(1 1)=((c56*x(12))+((c56+c45)*x(11))+(c45*x(10))+(K5*x(6))... 
~((K4+k5)*x(5))+(k4*x(4)))/J5; 

xdo(12)=(-Tload-(c56*x(12))+(c56*x(1 1))-(k5*x(6))+(K5*x(5)))/j6; 

xdot=xdo', % THIS IS THE VECTOR DEFINING THE EQUATIONS OF MOTION 


56 


LIST OF REFERENCES 


1. Wilson, W. K., “Practical Solution of Torsional Vibration Problems,” Third Edition 
Revised, pp. 499-661, John Wiley & Sons Inc., 1956. 


2. Kabele, D. F., “A New Approach in the Simulation of Crankshaft Torsional Vibration,” 
Society of Automotive Engineers Paper Number 844140, 1984. 


3. Citron S. J., O'Higgins J. E., Chen L. Y., “Cylinder by Cylinder Engine Pressure and 
Pressure Torque Waveform Determination Utilizing Speed Fluctuations,” Society of 
Automotive Engineers Paper Number 890486, 1989. 


4. Rizzom, G., “Diagnosis of Individual Cylinder Misfires by Signature Analysis of 
Crankshaft Speed Fluctuations,” Society of Automotive Engineers Paper Number 890884, 
1989. 


5. Sobel J. R., Jeremiasson J.,Wallin C., “Instantaneous Crankshaft Torque Measurement 
in Cars,” Society of Automotive Engineers Paper Number 960040, 1996. 


6. Brown T. S., Neill W. S., “Determination of Engine Cylinder Pressures from 
Crankshaft Speed Fluctuations,” Society of Automotive Engineers Paper Number 920463, 
noo 2. 


7. Bell, J. E., “Assessment of Diesel Engine Condition Using Time Resolved 
Measurements and Signal Processing,” Naval Postgraduate School, Monterey CA, M.Sc. 
Thesis, Sep. 1996. 


8. Detroit Diesel, Equivalent Mass Elastic System 3-53 Engine, Official Correspondence, 
May 1997, Detroit Diesel Corporation. 


9. Detroit Diesel Series 3-53 Crankshaft Dimensions, 1957, Detroit Diesel Corporation. 


10. Detroit Diesel Engines, Series 3-53 Service Manual, May 1990, Detroit Diesel 
Corporation. 


11. Detroit Diesel, Series 3-53 N. A. Performance Information and Peak Firing Pressure, 
Facsimile, Aug 1996, Detroit Diesel Corporation. 


12. Armstrong, R. A., “Fault Assessment of a Diesel Engine Using Vibration 


Measurements and Advanced Signal Processing,” Naval Postgraduate School, Monterey 
CA, M.Sc. Thesis, Dec. 1996. 


7 


13. SuperFlow, Inertia Data for Water Brake Dynamometer, Facsimile, Aug 1997, 
SuperFlow Corporation. 


58 


INITIAL DISTRIBUTION LIST 


Defense Technical Information Center 
8725 John J. Kingman Rd., STE 0944 
Ft. Belvoir, VA 22060-6218 


Py choy mio MET MrAgy eran este at ech cn scars vasvios<-cndSas s Segre vec euan vont ete iae ca tevis cesses. 2 
Naval Postgraduate School 
411 Dyer Rd. 

Monterey, Ca 93943-5101 


Wenanumentae malninanime Ode WWE ey. chy ce sshssex bss ie sess cs noe ees 
Department of Mechanical Engineering 

Naval Postgraduate School 

Monterey, CA 93943-5101 


Professor Knox T. Millsaps Jr., Code ME/MI........000...000 ecco ccccec cece ccceeeeeccceeeeneeeee 4 
Department of Mechanical Engineering 
Naval Postgraduate School 
Monterey, CA 93943-5101 


Oirrewlar Officer, Code 3450.20 oi.c. fe ete sees sacn ws oeoakicsl ee ee ee. 1 
Department of Mechanical Engineering 
Naval Postgraduate School 
Monterey, CA 93943-5101 


De molt, Wiesel CONPonariCm er eas ioas ie doce. s2esevesceevn evovapieearstesaees 1 
Mr. Kent Gregorius 

13400 Outer Drive, West 
Detroit, MI 48239-4001 


VANES COIN eA SQNED Oe ees os hasdiiewacislne ween cc écsheetedieccceadenevesevsovvesescsleanseemt 
880 Woodland Ave. 
Corydon, IN 47112 


aI WVU OM race hs. Shove cba dosocs obs cascduc cedcoelusaesedeesQlsesbuelsensesecntas 3 


313 Davis St. 
Portsmouth, RI 02871 


ae) 





DUDLEY KNOX LIBRARY 
NAVAL POSTGRADUATE SCH, 
MONTEREY CA 93943-5101 








WU 
32768 00342145 4 f 


—s 





